High Precision Measurement of the Thermal Exponent for the three-dimensional XY 

Universality Class 
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Simulations results are reported for critical point of the two-component cj> 4 field theory. The 
correlation length exponent is measured to high precision with the result v = 0.6717(3). This value 
is in agreement with recent simulation results [Campostrini et al, Phys. Rev. B 63, 214503 (2001)], 
and marginally agrees with the most recent space-based measurements of the superfluid transition 
in 4 He [Lipa et al, Phys. Rev. B 68, 174518 (2003)]. 
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Universality at critical points is among the most 
beautiful and powerful ideas to emerge from statistical 
physics yj. The unversality hypothesis asserts that crit- 
ical exponents and some other asymptotic properties of 
critical points are independent of microscopic details and 
depend on only a few properties such as the symme- 
try of the order parameter and spatial dimensionality. 
Universality permits the calculation of critical exponents 
for experimental systems using simplified and optimized 
model systems with the same symmetries. The theory of 
critical phenomena also predicts relationships amongst 
critical exponents. Among these relations are the hy- 
perscaling relation, a — 2 — dv, between the specific heat 
exponent a, the correlation length exponent v and the di- 
mensionality d and the Josephson relation, £ = (d — 2)i> 
for the superfluid stiffness exponent Q. Improvements in 
experiments, theory and computer simulations have led 
to increasingly strict tests of universality and scaling re- 
lations. By far the most accurate experimental measure- 
ments of critical exponents are for the superfluid tran- 
sition in 4 He, which is in the 0(2) or XY universality 
class. The purity of liquid Helium together with the sta- 
bility of temperature control and the accuracy of specific 
heat measurements at low temperatures means that the 
limiting factor in approaching the critical singularities is 
the rounding due to the Earth's gravitational field. To 
overcome gravitational rounding, space-based micrograv- 
ity experiments have been devised [3] that achieve four 
significant digit accuracy for the specific heat exponent a. 
These experimental results seemingly did not agree with 
either analytical or numerical calculations. The experi- 
mental value itself evolved with time due to reanalysis of 
the data, see Refs. jH3H The most recent reported 
value being v oxp = 0.6709(1), as inferred from the mea- 
sured value of a via the hyperscaling relation. 

A number of recent analytical and numerical calcula- 
tions were aimed at high precision determinations of v 
for the XY universality class. The vortex-loop calcula- 
tions by Williams [f| yield v — 0.6717. Torok and Hasen- 
busch 0] studied the two-component (j) 4 model via Monte 
Carlo simulations. This model has a parameter that ad- 
justs the softness of the potential constraining amplitude 
fluctuations of the order parameter. These authors took 



advantage of this freedom to suppress the leading correc- 
tions to scaling. This study resulted in v = 0.6723(3) (8) 
(the statistical and systematic errors are given in the first 
and second bracket, respectively). Later the effort was 
advanced by Campostrini et al. |8| who combined Monte 
Carlo simulations with a high-temperature expansion to 
obtain v = 0.67155(27), which agrees with the exper- 
imental result at the level of two standard deviations. 
Our study also employs a 4 model with fine-tuning of 
the Hamiltonian. 

Purely analytic studies of the XY critical point in- 
clude extensive treatments by Guida and Zinn-Justin, 
summarized in Ref. [jj, which yield v = 0.6703(15) 
for the perturbative seven-loop expansion in d = 3, 
and v — 0.6680(35) for the e-expansion up to the order 
e 5 . Jasch and Kleinert 01 developed a Borel resumma- 
tion technique in the context of a strong-coupling theory 
which results in v = 0.6704(7). The history of recent re- 
sults for the correlation length exponent is given in Fig.n 
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FIG. 1: Results for v as a function of time. Filled circles 
show the experimental values, open triangles depict the field- 
theoretic calculations, and squares show the Monte Carlo re- 
sults. 
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The purpose of this paper is to provide high precision 
simulations results for the exponent v for the 3D XY 
universality class using only the Monte Carlo simulations 
and to compare with other theoretical, computational 
and experimental results. 

We study a discrete 3D classical 4> 4 model, defined by 
the Hamiltonian 



H x 



(CV2)£>| 4 -mX> 



(1) 



where i and j label sites of a simple three-dimensional cu- 
bic lattice with periodic boundary conditions, (ij) stands 
for the pairs of nearest neighbor sites, and fa is the com- 
plex order parameter field. Since fa fields are continuous 
and unbounded, the model Q na s only two independent 
parameters: by rescaling fields one can always set one 
parameter to unity. Henceforth, we set fx = 1. 

In order to simulate the model Q , we employ the high- 
temperature expansion for the partition function which 
transforms the configuration space into that of closed ori- 
ented loops. The latter can be efficiently sampled by the 
worm algorithm pl| . which switches between the par- 
tition function and Green function sectors. The worm 
algorithm has virtually no critical slowing down and pro- 
vides a direct access to the statistics of winding number 
fluctuations, which, in turn, define the superfluid stiff- 
ness [T^: 



sL = J2(W*)/2d. 



(2) 



Here L is the linear system size, W a is the winding 
number in the a-th direction, and angular brackets de- 
note averaging over the Gibbs distribution. One can also 
devise direct Monte Carlo estimators for the derivatives 
of p s with respect to t and/or U, which involve cross- 
correlations between energy and winding numbers. 

In the renormalization group (RG) framework, the 
finite-size scaling of the superfluid stiffness obeys the re- 
lation 



p s L = f(x) + g UJ (x)y UJ L u + 



(3) 



where x = (L/^) 1 ^ is the dimensionless scaling variable, 
£ = £(t, U) is the correlation length, f(x) and g(x) are 
universal functions that are analytic as x — > 0, y w is the 
leading irrelevant scaling field, and dots represent the 
higher-order corrections. Field-theoretical calculations 
[J yield u = 0.802(11) (e-expansion) and u = 0.789(11) 
[d = 3 loop expansion), and the numerical analysis of 
Ref. @ gives w = 0.795(9). 

By differentiating Eq. © with respect to t for U = U c 
and then letting t — > t c , one transforms into 

R' c = AL 1/U {1 + CL~ U ) + higher-order terms. (4) 

Here the derivative of R = p s L is taken at the critical 
point, and A and C are non- universal constants. Eq. 



Q is especially convenient for the numerical data anal- 
ysis: the log-log plot of R' c versus L is a straight line 
with the slope \ jv for sufficiently large L. The second 
term adds a slight concavity or convexity to the curve 
for intermediate system sizes, depending on the sign of 
C. It is argued in Refs. [6,7] that there is little advan- 
tage in using improved models if corrections to scaling 
are included in the fits of MC data, and two alterna- 
tive strategies of dealing with the problem are suggested. 
One is that MC data are fitted by discarding correction- 
to-scaling terms, and the possible systematic error thus 
introduced is estimated from the universal amplitude ra- 
tios. The other is that MC data for the critical point are 
used in the analysis of the high-temperature expansion 
series. We demonstrate below that comparable accuracy 
can be achieved by Monte Carlo alone, using joint fits of 
several improved models. 

In order to locate a particular critical point, we fix 
some value of t and then plot R — p s L as a function of 
U for different system sizes. The crossing of these curves 
in the limit of L — > oo gives the critical point 0] . 

We first consider the critical point tj, — —0.0795548 
and U A = 0.4101562(14), which we denote by A. This 
critical point was previously studied by Campostrini et 
al. Q who performed an extensive and thorough search 
of improved models 01 • 

The accuracy of the critical point determination is ver- 
ified by an independent Monte Carlo measurement of the 
first and second derivatives of R. Indeed, suppose that 
for a given U — U c the value of t is off by St = t — t c . 
Then, expanding Eq. © up to 0(5t 2 ) and using the 
data from Table U one makes sure that (i) the devia- 
tions of R c from its universal value 0.2580(3) are con- 
sistent with St ~ 10~ 6 , and (ii) the 0(St 2 ) terms are 
smaller than statistical errors and, thus, can be safely 
neglected. For example, R A (L = 48) = 0.2583(1) and 
R'\(L = 48) = 2.65 x 10 4 . The derivatives for L ^ 48 
can be computed using R' oc L 1 ^ and R" oc L 2 ^, cf. 
Eq. ©. 

Table^lists the raw data for p' s L (see below for the dis- 
cussion of the dataset £>). Each data point was obtained 
from not less than 5 x 10 8 sweeps (upon equilibrating) 
and the accumulated data set represents approximately 
18 years of CPU time. Error bars are one standard devi- 
ation and were obtained using the blocking method. 

In order to magnify the fine details, we scale the nu- 
merical data with the experimental exponent, 



R'^L)^ 1 '^ 
R' c {L a )L~ 1/v - 



(5) 



and normalize it at an arbitrarily chosen value, Lq = 24. 
Figure |21 shows the data for critical point A, rescaled 
using Eq. JSJ. Given the error bars, it appears that 
corrections-to-scaling are relevant only for L < 10 and for 
L > 10 a straight line with the slope equal to \ jv — 1/ i/ exp 
is a good fit, see Eqs. 11151) . A linear fit of the data for 
L > 10 yields va — 0.67180(7), in a flagrant disagree- 
ment with the experimental value is cxp = 0.6709(1). Our 
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TABLE I: Results for critical points A and B. 



L 


R'a 


"•£5 


4 


2.0329(9) 


1.9906(3) 


5 


2.8414(5) 


2.7842(7) 


6 


3.7316(4) 


3.6583(9) 


7 


4.6955(5) 


4.6060(5) 


8 


5.7289(4) 


5.6215(9) 


9 


6.8265(7) 




10 


7.9848(9) 


7.839(1) 


11 


9.2031(13) 




12 


10.474(2) 


10.284(1) 


16 


16.074(4) 


15.784(2) 


20 


22.403(4) 


22.011(3) 


24 


29.396(7) 


28.883(3) 


32 


45.095(13) 


44.33(1) 


48 


82.48(3) 


81.05(2) 


64 




124.40(4) 


96 


231.56(17) 





value va is in agreement with the value v — 0.67155(27) 
of Rcf. 8] within the combined error bars. However, our 
new data for critical point A are more accurate and yield 
a smaller error bar. We have also studied a version of the 
improved link-current model |l5j |. which belongs to the 
same universality class. When these data were fit by a 
straight line, we observed even greater discrepancies with 




FIG. 2: The derivative of Ra rescaled via Eq. (JKJ. Dashed 
horizontal line corresponds to v = f oxp . See the text for 
discussion. 

To understand the apparent disagreement with ex- 
periment, we performed a series of simulations for crit- 
ical points in the vicinity of the critical point A to 
map out the region where C = C(t c , U c ) is small. We 
then carried out a large-scale simulation of one of the 
critical points, B, characterized by tg = —0.07142822 
and U B = 0.3605750(8). Table |T| shows the raw data 



for this critical point (R B {L = 48) = 0.2577(1) and 
R'£(L = 48) = 1.27 x 10 4 ). A naive linear fit of the 
dataset B for L > 10 yields v& = 0.67142(5), which is in- 
consistent with va- Critical points A and B are equally 
legitimate and we are left in a quandary as to whether to 
reject universality or to consider a more careful analysis 
of the data. We choose the later course. 

It becomes obvious that there is no easy way to im- 
prove the accuracy of the critical exponent calculations 
by ignoring corrections to scaling, even when they ap- 
pear very small. To reconcile the results for the datasets 
A and B one must include the subleading corrections in 
the fits. If this is done for each dataset separately, it 
leads to a large increase (almost by an order of magni- 
tude) in the uncertainty of the fits, as has been noted in 
Ref. @. Tighter error bars can be obtained performing 
a joint fit of two datasets. More specifically, we perform 
a 6-parameter fit according to 

lnq~B+(~-—}ln{L/L ) + CL-" ) (6) 

as follows from Eqs. (@J and © for C <C 1. The fitting 
parameter B is introduced in order to undo the artificial 
normalization q(L = Lq) = 1. Here each critical point 
has its own amplitudes B and C, while the universal ex- 
ponents v and us are shared between the datasets, which 
yields six fitting parameters. 

We performed the joint fit according to Eq. © via 
stochastic minimization of x-square. We constrained the 
exponent lu to within ±0.03 around 0.795, which is the 
3(7 interval of the established value 8] . We also discarded 
small system sizes L < L cuto s = 10. The optimization 
procedure yields v = 0.6717(3) at the confidence level 
[lj] CL « 0.43. This value is consistent with the result of 
Ref. [a, and also marginally agrees with the experimental 
value [3j . The fit is depicted in Fig.|3|and yields 

{C A = (1.5 ±0.5) x 10~ 3 , 
C B = (-1.0 ±0.3) x 10~ 2 , (7) 
v = 0.6717(3). 

Note that the dataset A points are exactly the same as in 
Fig- El and a slight curvature is visible in Fig. which is 
unaccounted for in Fig. [21 It is also worth noting that the 
best-fit value of u> is oj = 0.81(1) which is consistent with 
the field-theoretic estimates • The error bars in Eq. 
reflect both stochastic and systematic uncertainties, the 
latter being estimated via changing the £ C utoff for the 
dataset A and/or B. 

One might question the relevance of higher order cor- 
rections to scaling omitted from Eq. 0. The next order 
correction is proportional to the square of the leading ir- 
relevant field, oc C 2 L~ 2uJ , i.e. it contains the square of 
the already small amplitude C 0] . A sizable cumulative 
amplitude of higher-order corrections would exhibit itself 
for the smallest system sizes, L < 12. In the contrary, the 
scaling curves for both models under consideration have 
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FIG. 3: Data from Table Q rescaled via Eq. @. Circles: 
dataset A, diamonds: dataset B. Dashed horizontal line cor- 
responds to v = i^ CX p. Solid lines are the results of the joint 
fit of both datasets (see the text for explanation). Note the 
scale of the vertical axis. 



overall vertical scale on the order of 10 3 for system sizes 
as small as L = 4. 



The apparent controversy concerning the 3D XY uni- 
versality class is resolved. Our calculated value of v = 
0.6717(3) is in agreement with the results of recent sim- 
ulation |8( , and also marginally agrees with the best ex- 
perimental value i^cxp = 0.6709(1), deduced from a via 
hyperscaling relation. We demonstrate that one has to 
be careful when working with improved models and al- 
ways include corrections to scaling into the fit, even when 
the data set allows a good linear fit. The fourth-digit 
accuracy in critical exponents can be reached by simul- 
taneously fitting more than one critical point. 
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